package scu.maqiang.mesh.test;

import scu.maqiang.fes.BVPType;
import scu.maqiang.fes.FES2T31;
import scu.maqiang.mesh.Mesh2T3;
import scu.maqiang.mesh.ScalarFunc;
import scu.maqiang.numeric.Direct;
import scu.maqiang.numeric.NewIterSSolver;
import scu.maqiang.numeric.SRMatrix;

public class NACA0012 {

	public static void main(String[] args) {

		Mesh2T3 mesh = new Mesh2T3().ring(0.5, 1.0, 80, 20);
		mesh.translate(0.5, 0.0);
		FES2T31 fs = new FES2T31(mesh);
		SRMatrix A = new SRMatrix(fs.GetNdof());
		double[] RHS = new double[fs.GetNdof()];
		double[] NACA_Y = new double[fs.GetNdof()];
		fs.assembleStiff(new double[] {1.0}, BVPType.COMMON, A);
		fs.applyBC_MBN(A, Direct.All, 1, 2);
		ScalarFunc bdFuncOut = (xy, label, param) -> {
			double xx = xy[0];
			double y = xy[1];
			double curve1 = 0.17735 * Math.sqrt(xx) - 0.0755970 * xx - 0.2128360 * xx*xx + 0.173630 * xx*xx*xx - 0.062540 * xx*xx*xx*xx;
			return y > 0.00? curve1:-curve1;
		};
		ScalarFunc bdFuncInner = (xy, label, param) -> xy[1];
		ScalarFunc[] rf = {bdFuncOut};
		fs.applyBC_MBN(RHS, Direct.All, rf, null, 1);
		rf[0] = bdFuncInner;
		fs.applyBC_MBN(RHS, Direct.All, rf, null, 2);
		NewIterSSolver solver = new NewIterSSolver(A);
		solver.PCGSSOR(RHS, NACA_Y, 1.5e0, 1);
		//solver.GMRES_SSOR(RHS, NACA_Y, 1.5, 10, 2);
		for(int i = 0; i < mesh.getNv(); i++) {
//			mesh.nodes.get(i).coord[1] = NACA_Y[i];
		}
		mesh.toTecplot("NACA0012.dat");
		
	}

}





